module Auty
  implicit none
  type header_key
     integer*4 sizeof_hdr
     character*10 data_type
     character*18 db_name
     integer*4 extents
     integer*2 session_error
     character*1 regular
     integer*1 hkey_un0
  end type header_key
  type image_dimension
     integer*2 dim(8)
     character*4 vox_units
     character*8 cal_units
     integer*2 unused1
     integer*2 datatype
     integer*2 bitpix
     integer*2 dim_un0
     real*4 pixdim(8)
     real*4 vox_offset
     real*4 funused(3)
     real*4 cal_max
     real*4 cal_min
     real*4 compressed
     real*4 verified
     integer*4 glmax,glmin
  end type image_dimension
  type data_history
     character*80 descrip
     character*24 aux_file
     integer*1 orient
     character*10 originator
     character*10 generated
     character*10 scannum
     character*10 patient_id
     character*10 exp_date
     character*10 exp_time
     character*3 hist_un0
     integer*4 views
     integer*4 vols_added
     integer*4 start_field
     integer*4 field_skip
     integer*4 omax,omin
     integer*4 smax,smin
  end type data_history
  integer,parameter:: DT_NONE          =   0
  integer,parameter:: DT_UNKNOWN       =   0
  integer,parameter:: DT_BINARY        =   1
  integer,parameter:: DT_UNSIGNED_CHAR =   2
  integer,parameter:: DT_SIGNED_SHORT  =   4
  integer,parameter:: DT_SIGNED_INT    =   8
  integer,parameter:: DT_FLOAT         =  16
  integer,parameter:: DT_COMPLEX       =  32
  integer,parameter:: DT_DOUBLE        =  64
  integer,parameter:: DT_RGB           = 128
  integer,parameter:: DT_ALL           = 255     
  type ana_hdr
     type(header_key) hk
     type(image_dimension) id
     type(data_history) dh
  end type ana_hdr
contains

  integer function ana2mtv(fin,fout,tl)
    use suty
    ! fin = a.img, fout = a.mtv
    character*(*),intent(in)::fin,fout
    character*(*),intent(in)::tl ! toplabel
    character*512 fhdr,f
    character*32 cmax,cnx,cny
    integer lin,lout,istat
    type(ana_hdr) hdr
    integer*2,allocatable:: wbuf(:,:,:)
    integer*4,allocatable:: dbuf(:,:,:)
    real*4,allocatable::rbuf(:,:,:)
    real*8,allocatable::buf(:,:,:),b(:,:,:)
    integer nx,ny,nz
    integer ix,iy,iz

    ana2mtv=-1

    lin=len_trim(fin)
    fhdr=fin(1:lin-3)//"hdr"
    open(unit=1,file=fhdr,iostat=istat,status="old",access="stream")
    if(istat/=0) return
    read(1,iostat=istat) hdr
    close(1)
    if(istat/=0) return

    open(unit=1,file=fin,iostat=istat,status="old",access="stream")
    if(istat/=0) return

    nx=hdr%id%dim(2)
    ny=hdr%id%dim(3)
    nz=hdr%id%dim(4)
    allocate(buf(nx,ny,nz))
    allocate(b(nx+2,ny+2,nz))

    select case(hdr%id%datatype)
    case(DT_SIGNED_SHORT)
       allocate(wbuf(nx,ny,nz))
       read(1,iostat=istat) wbuf
       buf=wbuf
       deallocate(wbuf)
    case(DT_SIGNED_INT)
       allocate(dbuf(nx,ny,nz))
       read(1,iostat=istat) dbuf
       buf=dbuf
       deallocate(dbuf)
    case(DT_FLOAT)
       allocate(rbuf(nx,ny,nz))
       read(1,iostat=istat) rbuf
       buf=rbuf
       deallocate(rbuf)
    case(DT_DOUBLE)
       read(1,iostat=istat) buf
    case default
       istat=1
       return
    end select
    close(1)

    b=0
    do iy=1,ny
       do ix=1,nx
          b(ix+1,iy+1,:)=buf(ix,iy,:)
       end do
    end do
    deallocate(buf)

    if(istat/=0) return
    lout=len_trim(fout)
    cnx=i2a(nx+2)
    cny=i2a(ny+2)
    cmax=f2a(maxval(b))
    do iz=1,nz
       f=fout(1:lout-4)//trim(i2a(iz))//".mtv"
       open(unit=2,file=f,iostat=istat,status="replace",access="stream")
       if(istat/=0) return
       write(2) "$DATA=CONTOUR\n"
       write(2) "%toplabel="""//trim(tl)//"""\n"
       write(2) "%subtitle=Z="//trim(i2a(iz))//"\n"
       write(2) "%xlabel=""""\n"
       write(2) "%ylabel=""""\n"
       write(2) "%contstyle=2 nsteps=50\n"
       write(2) "%interpolate=3 bitmap=False\n"
       write(2) "%NX="//trim(cnx)//"\n"
       write(2) "%NY="//trim(cny)//"\n"
       write(2) "%XMIN=0 XMAX="//trim(cnx)//"\n"
       write(2) "%YMIN=0 YMAX="//trim(cny)//"\n"
       write(2) "%cmin=0 cmax="//trim(cmax)//"\n"
       write(2) "%BINARY\n"
       write(2) b(:,:,iz)
       write(2) "$END"
       close(2)
    end do
    deallocate(b)
    ana2mtv=0
  end function ana2mtv

end module Auty
!
